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Abstract 



The advantage of using a Discrete Variable Representation (DVR) is that 
the Hamiltonian of two interacting particles can be constructed in a very 
simple form. However the DVR Hamiltonian is approximate and, as a con- 
sequence, the results cannot be considered as variational ones. We will show 
that the variational character of the results can be restored by performing 
a reduced number of integrals. In practice, for a variational description of 
the lowest n bound states only n(n + l)/2 integrals are necessary whereas 
D{D + 1)/ 2 integrals are enough for the scattering states (D is the dimension 
of the S matrix). Applications of the method to the study of dimers of He, 
Ne and Ar, for both bound and scattering states, are presented. 



1 Introduction 



The solution of a quantum mechanical problem can be obtained using a 
variational principle. It provides second order estimates of binding energies 
and scattering matrices as well as the corresponding wave functions. For 
bound states the Rayleigh-Ritz variational principle can be used whereas for 
scattering states the Kohn variational principle (KVP), among others, can 
be employed. A common procedure is to combine the variational principle 
with the expansion of the wave function over a complete basis. Following 
Ref.pp we will refer to this method as the variational basis representation 
(VBR). The accuracy of such a method is strictly connected to the size of 
the basis set employed. In practical applications the basis set, which in 
many cases is infinite, is truncated. 

For the specific case of bound state calculations the implementation of 
the VBR method leads to the solution of a generalized eigenvalue problem. 
The corresponding eigenvalues and eigenvectors represent upper bounds to 
the exact energy levels of the Hamiltonian and first order approximations to 
the associated wave functions. If N is the dimension of the truncated basis 
the Hylleras-Undheim theorem guarantees that 

E X+1 (N) > E X (N) > E X (N + 1) , (1) 

and 

Jim E X (N) = e x , (2) 

where e x represents the exact eigenvalue of level A. As there is no approxi- 
mation involved in the calculation of the Hamiltonian matrix elements, the 
accuracy of the VBR is directly related to the completeness of the basis em- 
ployed. Accordingly, by increasing the dimension of the basis it is possible 
to obtain solutions extremely close to the exact ones. 

For scattering states the situation is slightly different. The wave function 
describing a scattering state is not an L 2 function, however its form outside 
the region where the collision takes place is in general known. Therefore the 
configuration space can be divided into two regions: the asymptotic region 
in which the particles are either free or interacting through a long range 
(coulombic) potential, and the internal region. In the asymptotic region the 
scattering wave function can be described in terms of Bessel or Coulomb 
functions. In the internal region the scattering state can be expanded over 
a complete L 2 -basis. The linear coefficients of the expansion and the S 
matrix elements giving the relative weights between the ingoing and outgoing 
asymptotic solutions can be obtained from the KVP. The implementation of 
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the Kohn variational principle using a VBR leads to a two-step procedure. 
The first step consists in obtaining a first order estimate of the scattering 
matrix S and of the wave function by solving a linear non-homogeneous 
system of equations of dimension N + D, with D the dimension of the S 
matrix. A second order estimate of the S matrix is obtained by inserting 
the first order solution in the functional: 

[S a f>] = S afi + i<*p\H-E\$+> , (3) 

with H and E the Hamiltonian and the energy of the system respectively. 
S a n is the first order estimate of the corresponding S matrix element and 
i s the outgoing (ingoing) solution corresponding to channel a (/?). 
The convergence of the second order estimate of the S matrix elements can 
be studied for increasing values of N. Although in this case there is no 
guarantee of a monotonous convergence, in contrast with the case for bound 
states, the properties of the KVP have been extensively studied in past 
years |2j. 

The key hypothesis in the variational theorems when implemented with 
the VBR method is that all the required matrix elements need to be cal- 
culated exactly, or at least to a very high accuracy. This means that the 
computational cost of the calculation increases as A 2 . It is quite common 
to choose orthogonal polynomials with appropriate weight functions as ba- 
sis sets, so that both the norm and the kinetic energy matrix elements are 
analytical, and the only numerical integrations required are those of the 
potential matrix elements. 

A different, but related technique, is the Discrete Variable Representa- 
tion (DVR) which is widely used in the description of molecular systems. 
Examples are the calculations of vibrational spectra, scattering problems or 
photo-dissociation processes [3J|I]. The DVR method is described in detail 
in a recent review 0. The basic property of the DVR is that local operators 
depending on the inter-particle distance (which include the potential) are 
diagonal in the DVR basis. However, in the case of the potential matrix 
the diagonal form is a result of an approximate treatment of the integrals 
defining the corresponding matrix elements. For example, in the case of an 
iV dimensional basis of orthogonal polynomials multiplied by appropriate 
weight functions, the approximation consists in replacing the integrals by 
the related N point Gaussian quadrature formula. In the limit of N going 
to infinity, both the VBR and the DVR approximations for the energy levels 
converge to the exact solutions. The advantage of the DVR with respect to 
the VBR is that there is no need to compute A 2 integrals, and the price to 
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pay is the loss of the variational character of the calculated eigenvalues and 
eigenvectors. 

In the present paper we would like to investigate the possibility of pro- 
ducing variational estimates using the DVR method. For bound states we 
proceed as follows. The solution of an N dimensional DVR problem yields 
N eigenvalues and the corresponding N eigenvectors -0^ VR . The low- 

est eigenvalue i^o >VR is an approximation to the exact ground state eo- Due 
to the absence of the variational character, i?Q )VR could be either greater 
or smaller than eo and this ambiguity remains even for large values of N. 
However the quantity 



is an upper bound to eo provided that the above integral is calculated accu- 
rately. To this aim, the wave function, which is known in the DVR 
basis, that is as amplitudes at the DVR points, has to be transformed to the 
r-space. As a result, only one numerical integral needs to be evaluated in or- 
der to obtain an upper bound to the ground state energy. As we will see this 
procedure can be extended to produce upper bounds to the lowest n levels 
belonging to a band of the exact Hamiltonian. In this case a minimum of 
n(n + l)/2 integrals have to be computed. In general the number N of basis 
functions required to describe correctly the first n states of an spectral band 
satisfies N » n, so that the number of integrals to be performed in order 
to produce the variational bounds is much smaller than those needed with 
the VBR. In cases where it is possible to obtain good solutions with N ~ n, 
the alternatives to the VBR method, such as the one here presented, do not 
have any particular advantage since the basis functions used to expand the 
wave functions are already close to the exact solutions. 

For scattering states a similar procedure can be followed. If 5 DVR is the 
first order estimate of the S matrix element calculated from the KVP using 
the DVR technique, and ■0 DVR is the corresponding scattering state, we can 
use eq. (J5J to calculate the second order estimate of the S matrix element. 
In order to produce a true variational estimate, the integral of eq. must 
be calculated accurately. This can be achieved again by transforming the 
DVR solution to the r-space. In this manner we have obtained a variational 
estimate of the S-matrix elements computing D(D + l)/2 numerical inte- 
grals, where D is the dimension of the S'-matrix. In most cases D is much 
smaller than the dimension of the basis N. 

In order to illustrate the application of the above formalism, we have 
studied dimers formed by two equal rare gas atoms. The spectrum of small 
clusters of He, Ne, Ar, Kr and Xe has been the subject of recent investiga- 




(4) 
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tions (see Ref. jS] and references therein). In particular, the study of the 
Helium dimer has proved very challenging, due to the extraordinary fea- 
tures of the He-He interaction that make the He2 molecule very elusive [7] . 
Difficulties in the theoretical study of rare gas clusters emerge from the char- 
acteristics of the atom-atom interaction, whose attractive part is of the van 
der Waals type and whose short range part consists of a hard repulsive core, 
causing those systems to be strongly correlated. As a result, the rotational- 
vibrational spectrum of these dimers differs substantially from the spectrum 
of typical covalent or polar molecules. The helium dimer has only one vibra- 
tional state, no rotational spectrum and its binding energy is about seven 
orders of magnitude smaller than that of traditional molecules. Some po- 
tential models predict three and nine vibrational states, respectively, for the 
Ne and Ar dimers however only some of those states have been observed. 
When the VBR is used to describe the bound states of these dimers a large 
number of basis functions is required to account for the strong atom-atom 
correlation. On the other hand the limited number of bound states makes 
the van der Waals dimers well suited for applications of the method outlined 
here. We have studied the dimers He2, Ne2 and Ar2 with state-of-the-art 
atom-atom potentials. In Section II we outline the method as applied to 
bound states. The convergence patterns of the energies of different states, 
as well as of other observables, are shown for comparison. In Section III we 
present results for the scattering lengths and the low energy phase-shifts. 
Section IV is devoted to the conclusions. 

2 Two particles bound states using VBR and DVR 

The center of mass Hamiltonian operator for two identical particles is written 



where M is the atomic mass and V(r) represents the interaction. In order to 
implement the VBR method we introduce the orthonormal basis $>kim( r ) = 
(frlk^Yimif) ■ The radial basis, which in the following we assume to belong 
to a family of orthogonal polynomials, satisfies 



For the sake of clarity, and without loss of generality, we can limit the 
discussion to I = states and call the corresponding radial basis elements 



as 



H = -—V 2 + V(r) 



(5) 




(6) 
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4>k (k = 0, 1, . . .). The radial part of the wave function corresponding to the 
level A is expanded in terms of the first N basis elements as 

N-l 

Mr)=Y, A kMr) ■ (7) 
fc=0 

The linear coefficients and the upper bounds E\ to the energy levels 
are obtained from the eigenvalue problem resulting from the Rayleigh-Ritz 
principle: 

N-l 

{H w ~ E x S w )Ay = (8) 

fe'=0 

Both the eigenvalues E\ and the corresponding eigenvectors Ay depend 
on the size N of the variational problem solved as is explicitly evident, for 
instance, from eq. lfl|). Therefore, a more complete notation would be E\{N) 
and Ah(N). However, in the following we will drop the index N to simplify 
the notation. 

The Hamiltonian matrix elements H^y are the sum of the kinetic and 
potential energy terms 




With an appropriate choice of the basis set, the integrals corresponding 
to the kinetic energy elements K^y can be obtained analytically whereas, in 
general, the integrals corresponding to the potential energy elements V^y are 
calculated numerically. As mentioned in the introduction, the accuracy of 
the numerical integration must be high enough to prevent error propagation 
in the computation of the eigenvalues of H. This meets the hypotheses of 
the variational theorem and assures convergence from above to the exact 
eigenvalues as the dimension N of the matrix increases. In the present work 
the relative accuracy in the computation of the potential energy elements is 
better than 1CP 7 . 

In an alternative approach, the DVR offers a computationally more ef- 
ficient method for obtaining good estimates of the eigenvalues and eigen- 
vectors of H . Here we briefly introduce the DVR for a basis of N orthog- 
onal polynomials times the appropriate weight functions JJ. In this case a 
one-to-one correspondence exists between the basis representation and the 



5 



representation in N Gaussian quadrature points, however the method is not 
limited to this case ^|EJ- The DVR corresponding to the set of N basis 
functions (f> k can be obtained from the following unitary matrix 

Tka = <t>k{x a )yfuj^ , (11) 

where {x a }, {u> a } (a = 0, . . . , N — 1) are the Gaussian points and weights 
corresponding to the quadrature formula 

N-l 

drr 2 4) k {r)V{r)4> k i(r) ph £ u a (fi k (x a )V(x a )(j) k >(x a ) . (12) 

The Gaussian points {x a } are the DVR points and they can be obtained 
as the eigenvalues of the N x N matrix representation of the coordinate 
operator r 

poo 

r kk > = / drr 2 cf> k (r)r 4> kl (r) . (13) 



o 

The corresponding eigenvectors form the T matrix defined in eq. (|11|) . 

The right hand side of eq. (|12|) defines the finite basis representation 
(FBR) approximation to V kk i and is called V^, BR . It corresponds to the 
computation of the potential matrix elements by quadratures. The potential 
energy operator in the DVR is defined as 

V™ = (TV FBR T) Q/3 = V{x a )5 a(] . (14) 

The DVR kinetic energy may be expressed as jf DVR = T l KT, where K is 
the matrix whose elements K kk i are given by eq. ©• Therefore we have 
introduced two isomorphic representations of H that have a very simple 
form. In _ff FBR = i^+y FBR ) the potential energy matrix has been calculated 
using an N point quadrature formula, whereas in _ff DVR = ^f DVR _|_ yDVR 
the potential energy matrix is diagonal. In both cases the kinetic energy 
matrix can be obtained analytically The two representations are related 
by the unitary transformation T, so they have the same set of eigenvalues 
j^FBRj = {^dvr^ T f ^fbr ig the e ig envector Q f H FBK with eigenvalue 

££ BR , then 

N-l 

^x BR (r) = E BX Mr) (15) 

k=0 

corresponds to the eigenvector of ff DVR in r-space, 

N-l 

^x BR (r) - ^ DVR (r) = £ C^\r) . (16) 

a=0 
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The coefficients in the above two equations are related to each other by the 
expression 

N-l 

Bl = J2 T ka C^ , (17) 

a=0 

moreover the DVR basis functions 0^ VR ( r ) in r-space are 

JV-l 

^° VR (r) = £ T kaMr) (a = 0,...,N-l) . (18) 

The r-space representation of the DVR basis element <^ VR (r) is a polyno- 
mial of degree N — 1 with zeros at the DVR points x@, with ft ^ a [jj. 

It is important to notice that the linear coefficients in eq. 1)15(1 differ 
from the VBR coefficients of eq. (JJJ) due to the quadrature formula used 
to calculate the potential energy matrix. As a consequence, the eigenvalues 
obtained cannot be considered as upper bounds. In fact, we will see that in 
several cases the energies £^ VR oscillate around the convergence value e,\. 
In order to obtain variational estimates, we apply the variational principle 
to the set n of DVR eigenvectors V , ? VR ( r )> i = 1 ; —,n, which approximate 
the n lowest levels of an spectral band of the exact Hamiltonian. In this 
paper we will discuss the case I = 0, however the method given below can 
be applied just as well to other /-bands. We can build the n x n matrix Ti 
as 

^■ = (^ DVR |^ + y|^ DVR ) , (19) 

computing the matrix elements of V numerically using the representation of 
^dvr j n r _ S p ace _ Thus, the computational effort of generating the matrix 
elements Hij is proportional to n(n + l)/2. Consequently, the eigenvalues 
£\ of Tiij represent real upper bounds to the eigenvalues e^. Moreover, as a 
consequence of the variational principle, the following relation holds: 

ex < E X (N) < £x{n) . (20) 

In fact, whereas E\(N) is the VBR eigenvalue calculated using N basis 
functions, £\(n) is the eigenvalue calculated using n specific combinations of 
the N basis functions and therefore it is obtained in a reduced Hilbert space. 
This procedure of generating variational bounds allows for a noticeable gain 
in computational time with respect to the VBR when the condition n « N 
is satisfied. In principle this condition could seem to be a restriction, but 
in many cases it is valid. For weakly bound molecules, as dimers of rare 
gases, the number of bound states is not very high. On the other hand, for 
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systems having a large number of bound states, the number of levels which 
are approximated well using N basis functions is generally smaller than N. 

In order to illustrate some applications of the above formalism we will 
calculate the bound states of dimers of He, Ne and Ar. For the He-He and 
Ne-Ne interactions we use the LM2M2 potential and the HFD-B potential, 
respectively, both proposed by Aziz and Slaman jlUl lllj . For the Ar-Ar 
system we use the HFD-C potential proposed by Aziz 12 j. The values 
h 2 /M = 43.281307, 8.584089 and 4.336093 Ka.u. 2 have been used for the 
He, Ne and Ar systems, respectively. With our choice of the inter-atomic 
potentials, and of the atomic mass values, the number of vibrational states 
that result is one for He2, three for Ne2, and nine for Ar2- The corresponding 
energy values are given in Table I. 

For the radial basis functions we have used the following orthonormal 
basis 

(2) 

where C k is generalized Laguerre polynomial depending on z = f3r with 
(3 a nonlinear parameter, which can be varied to improve the convergence 
patterns. 

The kinetic terms of eq. Q have been evaluated using the following 
analytical form: 

?2 * 2 r l„ -w, l,q 3,, 13. 



K - P h 

J\]chl 

M 



~{k + l){k + 2)8 kk , - ^k 3 --k 2 -^k-l 



(k < k'). 



(22) 

The non-linear parameter (5 was chosen to be 1 a.u. , 5 a.u. 1 and 10 
a.u. -1 for He2, Ne2 and Ar2, respectively. 

In figure 1 the convergence of the ground state energy of He2 is shown 
as a function of N for the three methods, VBR, DVR and the mean value 
calculated using the DVR wave function which in the following has been 
called (DVR). It may be observed that the DVR energy E7g R (dashed line) 
oscillates around the VBR energy Eq (solid line) even at very high values 
of N. The solution becomes stable for N > 500 revealing the particular 
structure of the He2- On the other hand the (DVR) energy £q (dotted 
line) has a much more stable pattern of convergence and already at N = 
200 its value coincides to five digits with the VBR result. In fact, in the 
figure, the VBR and (DVR) curves are almost indistinguishable. The (DVR) 
calculations can be used to estimate the quality of the DVR wave function. 
In Table II the total energy, the potential energy and the square root radius 
are given as functions of N for the three methods. As already mentioned, 
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we observe a better stability in the VBR and (DVR) energies. Conversely 
the DVR and the (DVR) potential energy and square root radius show small 
differences with respect to the VBR values. 

In figure 2 the patterns of convergence of the energy are shown for four 
selected vibrational states of A12 and Ne2, as a function of the number ./V 
of DVR points. For Ne2 the HFD-B potential predicts three bound states, 
however the highest excited state is very loosely bound and has been initially 
excluded from the analysis. Therefore the wave functions V? VR (A = 0, 1) 
have been used to calculate the matrix elements of the Hamiltonian as de- 
fined in eq. (|19|) . Accordingly an eigenvalue problem of dimension n = 2 
has been solved. For Ar2 the HFD-C potential predicts nine bound states, 
but again the highest excited state has been initially excluded, therefore an 
eigenvalue problem of dimension n = 8 has been solved. We have obtained 
a good convergence for all the state analyzed for both molecules, but for 
clarity here we limit the figure to show the energies corresponding to the 
ground states of the two molecules, the first excited state of Ne2, and the 
eighth level of Ar2. The convergence pattern of these states present some 
characteristics that merit discussion. As expected the VBR has a very stable 
convergence in all cases. The DVR energies approach the VBR result with 
a very fast convergence yielding oscillations of decreasing amplitude as N 
increases. Though the DVR results are qualitatively similar in the four cases 
presented, and also in all the other states examined, the (DVR) results for 
the two ground states have different patterns. Whereas £q for Ne2 (dotted 
line in panel (a)) presents a marked oscillatory pattern, £q for Ar2 (dotted 
line in panel (c)) presents a convergence extremely close to that of the VBR 
case. Furthermore, the two (DVR) excited states shown in figure 2 as dot- 
ted lines in panels (b) and (d) also present a marked oscillatory pattern. 
This non monotonic convergence is not in contradiction with the variational 
principle since the relation presented in eq. (j2()|) holds for each value of N. 
Moreover these oscillations are a consequence of the non-variational charac- 
ter of the DVR wave functions since in many cases we found that the N + 1 
result is worse than the iV result. 

It is possible to improve the convergence of the (DVR) by increasing the 
number of DVR functions included in the variational problem of eq. (|19|) . If 
all the N DVR wave functions are included in eq. (|19|) . the VBR and (DVR) 
methods become isomorphic. Whereas a minimum of n basis functions is 
sufficient to produce upper bounds to the n lowest states of a band, it is 
reasonable to expect that on increasing the number of functions the 
convergence pattern will improve. Let us call this number n' . In figure 3 
the convergence of the two energy levels of Ne2 (panels (a) and (b)) and 
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the eighth level of Ar2 (panel (c)) are given as a function of the number N 
of DVR points, for different sizes n' of the (DVR) problem. We consider 
the cases n' = 4,8,12 for Ne2 and n' = 12,20,30 for Ar2, corresponding 
to the solid, dashed and dotted-dashed lines, respectively. In all cases the 
plot of the VBR results is shown as a dotted line. The eighth level of Ar2 
represents the most difficult case for the (DVR) method as it requires n' 
to be consistently larger than n in order to have a smooth convergence 
pattern. However, we observe that there is a significant improvement in 
the convergence pattern with choices of n' that still satisfy the important 
condition n' <C N. Moreover, as expected from the inequality of eq. 1)20(1 . 
the E\ energies are always greater than the correspondening E\ energies. 

3 scattering states of two particles using VBR and 
DVR 

In order to produce variational estimates for the scattering matrix using the 
DVR method we use the Kohn variational principle in its general form [2]. 
Here we give a brief introduction to the method, limiting it to the case of 
local central potentials. Without loss of generality, we treat specifically the 
case I = and the case in which the collision proceeds along one open chan- 
nel. In this case the scattering matrix is a scalar quantity. The formalism 
can be easily generalized to treat more channels and the case in which the 
long range Coulomb potential is present |13| . 

The radial scattering wave function corresponding to a process at energy 
E can be written as a sum of two terms: 

^+(r)=^ c (r)+^+(r) . (23) 

The first term, ip c is the internal part and describes the system when the 
two particles are close to each other. It can be expanded in terms of N L 2 
basis functions just as for the bound states 

jV-l 

Mr) = E A k<t>k(r) ■ (24) 

k=0 

The second term describes a general asymptotic scattering state of two par- 
ticles and is defined as 

i;+(r)=n (r)+£n 1 (r) , (25) 
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where 



= y -^[uooMqr) +u 01 y (qr)] , (26) 
I Mq 

fti = y ^Koiown + uumqr)] . (27) 

are asymptotic scattering states, with q 2 = j^E. They are given as com- 
binations of the regular spherical Bessel function jo and the product of the 
irregular Bessel function yo an d a regularizing factor, namely 

y (qr) = (1 - e^ r )y Q {qr) . (28) 

The specific form of the regularizing factor is not crucial provided that the 
regularization is made in the internal region and yo — ► yo outside the range 
of the interaction. The quantity C expresses the relative weight of the two 
scattering asymptotic states Qq and f2i. The coefficients Uij form a matrix 
that can be chosen in accordance with the different meanings of the quantity 
C. For example the choices 



(29) 
(30) 



define the reactance matrix R and the scattering matrix S, respectively. 
The generalized KVP states that the functional 

[ £ i = £ -d^~'^' +> 1311 

is stationary with respect to variations of the parameters used to construct 
the wave function. We define H to be H = (H — E), det(n) is the de- 
terminant of the matrix u and ip~ is the complex conjugate of vp + . The 
normalization of the asymptotic states is defined so as to satisfy 

(Ool-HM) - (^i\H\n ) = ^det(-u) . (32) 




The unknowns in the wave function of eq. (|23|) are the N linear coeffi- 
cients Ak and the quantity C. The variation of the Kohn functional with 
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respect to these unknowns leads to the following N + 1 inhomogeneous sys- 
tem of equations 



N-l 

£ A k (<j> k \H\<f> k ,) -£(0 fc ,|H|Oa) = -(cf> k ,\H\no) , (k' = 0, . . . ,N — 1 

k=0 
N-l 



ZM*>\aw-afk\ttw = - det( " ) + 2( "' l T° >+2( "° |fl > > 

fe=0 



The coefficients A k as well as the first order estimate of the scattering matrix 
C are obtained from the solution of the above system of equations. The 
second order estimate [C] is calculated substituting the first order solution 
in eq. (|31j) and gives 

[C] = £--^r T f(Oo|y|fio) + £ 2 (Oi|^|^i) + £[(^o|i?|Oi) + (Oi|^|^o)] 
det(u) i- 

N-l N-l 

+ 2CY,M<t>k\H\Sli) + 2Y,A k (4> k \H\Sl Q )) . (34) 

k=0 k=0 

The convergence properties of the KVP variational principle have been 
extensively studied ^JE]- Occasionally singularities occur in its solution, 
however the complex form (£ = S) has a much more stable convergence 
pattern. Here we use both forms and check the consistency of the results by 
means of the relation 

S = (l + ifyil-iR)- 1 (35) 

which holds only for the exact matrices. Moreover the unitary condition for 
the S matrix, SS' = I is checked as an indication of the completeness of 
the basis used for the expansion of the internal part of the wavefunction. 

The DVR method can be used to calculate the first order estimate of 
the scattering matrix. For this we write the linear system of eq. (|33j) in the 
compact form 

HA = h , (36) 

where H is the (N + 1) x (N + 1) matrix defined to be H kk / = {(j)k\H\(j)k>) 
(k,k' = . . . iV) and writing (fix = Q\. Accordingly, the last element of the 
vector of coefficients A is An = — C and the elements of the vector b are 
defined to be 

b k = -(<f, k \H\Qo) k < N , (37) 

-det(u) + 2(fii|#|fio>+2(fi |#|ni) 
b N = -. • (38) 
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The latter can be reduced to bjy = (Qq\H\Qi) utilizing eq. (|3*21 

We can now extend the unitary transformation T which acts in the NxN 
space to the (N + 1) x (N + 1) space by defining the following unitary matrix 
f 

f ij = Tij i, j = 0. . .N- 1 , 

T Nj=f jN = j = . . . N - 1 , (39) 

T NN = 1 • 

As in the case of the bound state, we first introduce H FRR and b FBR in 
which the integrals involving the potential energy operator are calculated 
using an N point quadrature formula. The integrals involving the asymp- 
totic functions Q p are also calculated using quadratures. The corresponding 
vector of linear coefficients is A FBR . The DVR representation is defined to 
be 

^DVR A DVR = b DVR ^ (4Q) 

with 

H DVR = f l H FBR f , (41) 

A DVR = ^ A FBR (42) 
b DVR = ^ b FBR (43) 

The structure of the symmetric matrix i^ DVR is the following 

Hf™ = ^DVR + v{x . )6ij h j = N _ x f (44) 

H%] R = ^w][HQx] Xj j = 0...N-l, (45) 

JV-l 

HnT = E wMx a )[HQi] Xa ■ (46) 

a 

The structure of the inhomogeneous term b DVR is 

bf VR = -^U-[HQ ] Xj j = 0...N-l (47) 
JV-l 

^ VR = - E w a n (x Q )[Hni] X( , , (48) 

a 

where [HQ p ] Xa indicates the evaluation of the function HVL p at the DVR 
point x a . Prom the above relations we see that only two quadratures have 
to be performed in order to construct the DVR system of equations. The 
computational effort is limited as compared to the VBR technique. 
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The solution £ DVR = —^dvr re p resen ts the DVR first order estimate 
of the scattering matrix. It coincides with the solution £ FBR = _^4FB R ob- 
tained by solving the FBR system. However both differ from the first order 
scattering £, obtained by solving the system of eqs. using the VBR, 

due to the quadratures introduced to calculate the matrix elements. It is 
important to notice that the asymptotic integrals and (£Iq\H\Qi) , 

involved in the construction of the matrix element H^J^ and in the inhomo- 
geneous term 6^ VR respectively, do not depend on the radial basis 4>k used 
to expand the scattering wave function. When the VBR technique is used, 
these integrals are calculated numerically using standart techniques. In the 
DVR, as it is clear from eqs. (|46j) and (|48ft . they depend on the dimension N 
of the basis. In some cases, the use of the VBR values in H^ffi an< ^ kjv VR i 
instead of those obtained from the quadrature formula, fails to yield the 
solution (an example will be given at the end of Section III-A). Conversely, 
the use of the formulas of eqs. (|4*6|) and (fl8|) produces a convergence pattern 
similar to that of the bound states. This means that one has to be very 
careful in mixing exact integrals and quadratures when using the DVR. 

In order to produce a variational estimate of the scattering matrix we 
construct the DVR radial wave function in r-space, 

^DVR (r) = ^DVR (r) + fio ( r) + £DVR 0l ( r) ^ ^ 
where the internal wave function is 

N-l 

^(r) = £ ^ VR ^ VR (0 • (50) 

Using eq. (|BT)|) and eq. (|HTj) it is possible to calculate the second order 
estimate of the scattering matrix: 

[£] DVR = £ DVR_ 2^ DVR ^ DVR 

det(u) 

The integral A = (?/> DVR |.ff |-0 DVR ) converges since if |-(/> DVR ) — > for r > rj, 
where 77 is the interaction range. 

We now consider [£] DVR to be a second order variational estimate of 
the scattering matrix. It has been obtained calculating one integral, namely 
the integral A. For the case in which D channels are open, D(D + l)/2 
integrals must be evaluated. Since, in general, D <C N the computational 
effort needed to produce the DVR variational estimate is much smaller than 
that of the corresponding VBR. In the following subsections we discuss the 
results for zero and positive energy scattering separately. 
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3.1 zero energy case 

For zero energy scattering the asymptotic wave functions and the KVP 
assume a particular form. The zero energy wave function is 

Mr) = Mr) + Mr) , (52) 

with 

N-l 

Mr) = E A kMr) (53) 
k=0 

Mr) = fio(r)-afii(r) • (54) 

The asymptotic functions are defined to be S7o = <J^i an d ^i = 

V Zh?(^~ ~ e_7r )/ r an( i a 1S the scattering length. The Kohn functional for 
the scattering length can be obtained from eq. (|31jl in which C = R and 
taking the limit 

lim &U = —a (55) 
g^o q 

where 5 is the phase shift and we have used the definition R = tan 5. We 
obtain 

[a] =a + 2(MH\M ■ (56) 

The minimization of the above functional with respect to the set of coeffi- 
cients {Ak} and the scattering length a leads to a linear system of equations 
formally equal to that of eq. (J36|) or, in the case of the DVR, to that of eq. 
(|40j) . We solve both systems to calculate the VBR and DVR first order scat- 
tering lengths for the He2, Ne2 and Ar2 systems. As mentioned before, the 
second order estimates are obtained by substituting the first order solution 
in eq. (|5li|) calculating the integral (tpo\H\ipo) numerically. Again we use 
the label (DVR) for the second order estimate when using the DVR wave 
function. 

In figure 4 (left panel) we show the He-He scattering length (second or- 
der) calculated using VBR (solid line) and (DVR) (dashed line) as a function 
of the number N of basis functions. The converged quantity is a = 189.518 
a.u.. The attempt to correct the first order DVR result by calculating the 
second order estimate using quadratures produces a non-convergent value 
as is shown in the figure with the dotted line. In the following we call this 
result DVR g . To further analyze the differences between the first and second 
order, in the right panel of figure 4 we show the VBR and DVR first or- 
der results. Whereas the first order VBR result converges as the dimension 
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of the problem increases, the first order DVR result oscillates around the 
exact value even for N > 300. By comparing the first order DVR result, 
the DVRq result and the (DVR) result, we conclude that the error intro- 
duced by the quadratures in the integral (^q|J3"|^q) is of the same order of 
magnitude as the integral itself, therefore it cannot give the correction prop- 
erly. Conversely, using the DVR as a trial wave function the application of 
the Kohn functional naturally produces a better approximation. Finally, in 
figure 5 we show the convergence as a function of N of the (second order) 
scattering lengths for the Ne-Ne system (a) and Ar-Ar system (b). The 
converged values are a = 28.411 a.u. and 520.22 a.u., respectively. The 
solid, dashed and dotted lines represent the VBR, (DVR) and the DVRq 
second order estimates, respectively. For Ne the three calculations converge 
reasonably well. On the contrary, in the case of Ar the VBR presents a fast 
convergence, the (DVR) converges after marked oscillations and the DVR g 
fails to converge even for N ~ 200. This is a good example of where the 
corrective term proportional to A cannot be calculated using quadratures. 
The DVR wave functions for the He-He, Ne-Ne and Ar-Ar systems are given 
in figure 6 corresponding to cases N = 300, N = 200 and N = 400 respec- 
tively. Each function presents a number of nodes equal to the number of 
bound states supported. For the Ar-Ar system a basis of big dimension is 
needed (N = 400) in order to produce a wave function orthogonal to the 
very loosely E% state. The values used for the non-linear parameters are 
(3 = 1 a.u. -1 and 7 = 0.1 a.u. -1 for He, (3 = 5 a.u." 1 and 7 = 0.1 a.u. -1 for 
Ne, and (3 = 5 a.u. -1 and 7 = 0.5 a.u. -1 for Ar. 

Finally we would like to discuss the definition of the elements H^ffi 
and 6° VR given in eqs. (|46[) and (|48|) . They correspond to the asymptotic 
integrals (£li\H\£li) and (Slol-^l^i) obtained using quadratures and, because 
of that, they vary with N. Therefore at each successive step the DVR problem 
has been solved using different values for those elements. On the other 
hand, using the VBR, these elements has been computed numerically and 
independently of N. The differences can be small but significant, for example 
in some of the present calculations the relative difference between both types 
of integration is about 10 -5 , for N ~ 200. Analyzing the scattering length 
convergence pattern as a function of N using both types of integration we 
observed small differences in the cases of He and Ar. However the quadrature 
formula is to be preferred since it provides a first order estimate slightly 
closer to the final result. In the case of Ne the use of the quadrature formula 
is mandatory, since otherwise the first order estimates are far from any 
reasonable value, even for large values of N. This may be understood by 
comparing the scattering lengths with the interaction range 77 (estimates 
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of n are given below in B). For He and Ar the scattering length a is much 
larger than the interaction range whereas this is not the case for Ne. In 
the latter case there is a delicate cancellation between the internal and the 
asymptotic part of the wave function in the region 77 < r < a. Using DVR 
this cancellation can be properly achieved only when the definitions of eqs. 
£HJ) and (@SJ) are used. 



3.2 positive energy case 

Here we study the convergence properties of the I = phase shift 5 at very 
low energies using the DVR. In this energy region the phase shift and the 
scattering length are related by the effective range expansion 

+ ~7rQ > ( 57 ) 



tan 5 a 2 
where tq is the effective range parameter defined to be 



9 i-oo o*2 

r = 4 / [(r-a) 2 -—(r^ f]dr , (58) 
ar Jn M 



with ^0 given in eq. (|52j) . Since (ripo) — > J Mz(r — a) for values of r > rj, 



2% i 

the above integral goes very rapidly to zero outside the interaction region. 
Therefore ro is itself a measure of the interaction range rj. We obtained 
the values of 13.94 a.u., 15.65 a.u. and 59.25 a.u. for He, Ne and Ar, 
respectively. 

In Table III we study the convergence of the phase shifts 5 for the He-He 
and Ne-Ne systems at 50 mK and for the Ar-Ar system at 2 mK. The VBR, 
DVRg and (DVR) results correspond to second order estimates obtained 



replacing the solutions of the linear systems of eqs. (|36l40p in the Kohn 
functional of eq. (|31[) . Here we do not want to discuss the convergence 
properties of the KVP, which can be found in Refs. |14l I15j . However, in 
order to avoid possible singular solutions, we have applied the KVP to the 
two cases, C = R and £ = S, and checked the equivalence of the results in 
accordance with eq. (|35|) . In addition we have verified that on increasing 
the number of basis functions the quantity \SS^ — I\ ~ 10~ 8 in all the 
cases considered, which is of the same order of magnitude as the differences 
between the results using C = R or C = S. From the table we can see that 
the DVRq phase shifts have a slower convergence pattern in all cases. For 
the He-He system we have seen that at zero energy the DVR g result is not 
stable even for large N values. At 50 mK the situation is slightly improved 
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though the pattern of convergence of DVR g is still not satisfactory. For 
Ne and Ar we observe that big values of N are necessary to obtain stable 
results. This is due to the presence of excited states close to threshold which 
are generally difficult to describe, as in the case of the Eg state of Ar. 

In figure 7 we compare the quantity q/tan.5 (filled circles) to the r.h.s 
of eq. (|57|) calculated using the zero energy results (solid line). From this 
plot it is possible to extract the energy at which the phase shift starts to 
deviate from a linear relation. Below this energy the collision is not very 
sensitive to the details of the potential since the dynamics is governed by two 
parameters, the scattering length and the effective range. Greater sensitivity 
to the potential appears above this energy. For the He system the deviation 
from linearity starts at about i? = 0.1K whereas for Ne it appears at about 
E = 0.05 K and for Ar around E = 0.002 K. 

4 Conclusions 

In the present paper we have studied the possibility of using DVR wave 
functions to produce variational estimates to binding energies and scatter- 
ing matrices. The main advantage in using the DVR in the description of a 
quantum mechanical problem is the simplicity in constructing the Hamilto- 
nian matrix. However the DVR eigenvalues do not represent upper bounds 
to the exact levels. Therefore, in order to obtain estimates to the levels, it 
is necessary to produce a convergence pattern in terms of the number N 
of DVR points and analyze the stability of the results as N increases. In 
the case of the He-He system we have seen that the DVR eigenvalue oscil- 
lates around the exact level even for high values of N. On the contrary, the 
variational estimate obtained with the DVR wave function shows a better 
stability. As a result, calculating only one integral namely the mean value 
of the Hamiltonian with the DVR ground state wave function, we have im- 
proved the prediction using the DVR method. 

In order to extend the discussion to excited states, we have studied 
the Ne-Ne and Ar-Ar systems. We have shown that we need to perform 
at least n(n + l)/2 integrals to produce variational estimates of the first n 
levels. Although the variational estimates are always upper bounds, we have 
noticed that, for some levels, they do not converge monotonically, presenting 
instead oscillations. In order to improve this behavior we have considered a 
variational problem with an enlarged number (n') of DVR wave functions. 
We have seen that with n' ~ 2n we obtain convergence patterns very close to 
those obtained with the VBR. In the case in which the description includes 
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many excited states with a big spread in energy, we have found it necessary 
to further increase n'. Both n and n' are in general much smaller than N, so 
the possibility of generating variational bounds using the DVR is convenient 
as compared to the VBR, from the point of view of the computational effort 
required. 

For scattering states the Kohn Variational Principle gives a natural con- 
text in which the DVR wave function may be used to calculate second order 
variational estimates of the scattering matrix. The case of the He-He sys- 
tem at zero energy is a good example. The first order scattering length 
calculated using DVR oscillates around the exact value. The second order 
variational estimate calculated using the DVR wave function has a very fast 
convergence. Conversely, the second order when calculated using quadra- 
tures fails to reproduce the correct value for N < 500. This is a consequence 
of the particular structure of the He-He interaction. For the Ne-Ne system 
the differences between the second order calculated variationally or by using 
quadratures is not so pronounced as N increases. At the same time, for the 
Ar-Ar system the second order when calculated using quadratures fails to 
converge even for values of N of the order of 400. Perhaps this is the best 
example in which the application of the KVP to correct the first order es- 
timate is clearly necessary. In general we have observed that the scattering 
calculations need a bigger basis set to obtain converged results. This follow 
from the fact that the Ne-Ne and the Ar-Ar interactions predict an excited 
state close to zero energy. 

Finally we would like to mention the possibility of extending the present 
study to the three-body problem. In this case the dimension of the matrices 
needed for the description of the processes are much bigger. So the capability 
of producing convergence patterns in terms of N could be limited. Hence 
variational estimates such as those presented here could help to improve the 
DVR predictions. A study of this subject is at present in progress. 
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Table 1: Spectrum of the He, Ne and Ar dimers calculated with the LM2M2, 
HFD-B and HFD-C potentials, respectively. 





He(mK) 


Ne(K) 


Ar(K) 


Eq 


-1.3020 


-24.4422 


-121.8571 


Ei 




-4.5281 


-84.9229 






-0.0327 


-55.4794 


E% 






-33.1574 


E 4 






-17.4289 


E 5 






-7.5579 


Eq 






-2.3509 


E 7 






-0.3409 


Eg 






-1.9 10" 4 



Table 2: The He2 binding energy, potential energy and mean square root 
radius, calculated with the VBR, DVR and (DVR) methods, as a function 
of the size N of the basis set employed. Energies are in mK, the radii in a.u. 





VBR 


DVR 


(DVR) 


N 




E 


V 


vV> 




E 


V 


vV> 




E 


V 


vV> 


50 


-0 


0588 


-141.786 


67.34 


-0 


0742 


-140.622 


67.65 





1002 


-140.448 


67.65 


100 


-1 


2261 


-108.303 


106.66 


-1 


2337 


-108.480 


106.64 


-1 


2177 


-108.464 


106.64 


150 


-1 


2943 


-102.049 


125.33 


-1 


2938 


-102.034 


125.34 


-1 


2942 


-102.035 


125.34 


200 


-1 


3012 


-100.904 


131.86 


-1 


3026 


-100.958 


131.80 


-1 


3012 


-100.957 


131.80 


250 


-1 


3019 


-100.718 


133.62 


-1 


3025 


-100.742 


133.59 


-1 


3019 


-100.742 


133.59 


300 


-1 


3020 


-100.690 


134.01 


-1 


3009 


-100.648 


134.06 


-1 


3020 


-100.649 


134.06 


350 


-1 


3020 


-100.686 


134.09 


-1 


3027 


-100.714 


134.05 


-1 


3020 


-100.713 


134.05 


400 


-1 


3020 


-100.686 


134.10 


-1 


3014 


-100.662 


134.13 


-1 


3020 


-100.662 


134.13 


450 


-1 


3020 


-100.686 


134.10 


-1 


3015 


-100.665 


134.13 


-1 


3020 


-100.665 


134.13 


500 


-1 


3020 


-100.686 


134.11 


-1 


3019 


-100.682 


134.11 


-1 


3020 


-100.682 


134.11 
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Table 3: Convergence of the second order estimates of the phase shift 5 (in 
radians) at different energies for the He-He, Ne-Ne and Ar-Ar systems, as a 
function of the dimension iV of the basis. The calculations using the VBR, 
the DVR g and the (DVR) methods are compared. 



N 


VBR 


DVR 


(DVR) 




He-He (E 


=50 mK) 




20 


1.2563 - 


1.3694 ■ 


-1.5563 


40 


1.4685 


1.4476 


1.4627 


60 


1.4860 


1.4892 


1.4741 


80 


1.4767 


1.4635 


1.4784 


100 


1.4883 


1.4874 


1.4877 


120 


1.4883 


1.4884 


1.4882 


140 


1.4883 


1.4879 


1.4883 


160 


1.4883 


1.4884 


1.4883 


180 


1.4883 


1.4882 


1.4883 


200 


1.4883 


1.4884 


1.4883 




Ne-Ne (E 


=50 mK) 




20 


1.0573 


1.5592 • 


-1.5708 


40 


1.3272 


1.3450 - 


-0.0693 


60 


1.3309 


1.3456 


1.3209 


80 


1.3310 


1.3311 


1.3309 


100 


1.3310 


1.3314 


1.3310 


120 


1.3310 


1.3310 


1.3310 


140 


1.3310 


1.3310 


1.3310 


Ar-Ar (E=2 mK) 


10 


1.2477 - 


0.1641 


1.4324 


50 


-0.4193 - 


0.6421 - 


-1.5435 


100 


0.9310 


1.1382 


-1.4981 


150 


1.0099 


1.0077 


0.9992 


200 


1.0109 


1.0099 


1.0108 


250 


1.0110 


1.0107 


1.0111 


300 


1.0110 


1.0108 


1.0111 


350 


1.0110 


1.0108 


1.0111 
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Figure Captions 
Figure 1. 

Ground state energy of He2, calculated using the VBR (solid line), 
DVR(dashed line) and (DVR) (dotted line) methods, as a function of the 
dimension ./V of the basis. 

Figure 2. (a) ground state energy of Ne2, (b) first excited state energy of 
Ne2, (c) ground state energy of Ar2 and (d) eighth level of Ar2, calculated 
using VBR (solid line), DVR (dashed line) and (DVR) (dotted line), as 
a function of the dimension N of the basis. For clarity the results have 
been drawn with a solid line (VBR), a dashed line (DVR) and a dotted line 
((DVR)). 

Figure 3. Convergence patterns for the energy of the ground state of Ne2 
(a), the first excited vibrational state of Ne2 (b), and the eighth level of Ar2 
(c) , for different choices of the size n' of the restricted variational problem of 
eq. (fT9]l. For Ne2 the lines correspond to n' = 4 (solid), n' = 8 (dashed) and 
n 1 = 12 (dotted-dashed). For Ar2 the lines correspond to n' = 12 (solid), 
n' = 20 (dashed) and n' = 30 (dotted-dashed). The VBR results are shown 
as a dotted line for reference. 

Figure 4. Second order estimate for the He-He scattering length a calculated 
with the VBR (solid line), (DVR) (dashed line) and DVR g (dotted line) as 
a function of N (left panel). First order estimate for a, as obtained with the 
VBR method (solid line), and DVR (dotted line) as a function of N (right 
panel). 

Figure 5. Second order estimate for the Ne-Ne scattering length a calculated 
with the VBR (solid line), (DVR) (dashed line) and DVR q (dotted line) as 
a function of N (left panel). Second order estimate for the Ar-Ar scattering 
length a calculated with the VBR (solid line), (DVR) (dashed line) and 
DVRg (dotted line) as a function of N (right panel). 

Figure 6. Zero energy DVR wave functions for He2 (upper panel), Ne2 
(middle panel) and Ar2 (lower panel) calculated using ./V = 300, N = 200 
and N = 400 respectively. 

Figure 7. The r.h.s. of eq. (|57|) (solid line) compared to q/t&nS (filled 
circles) for He- He (upper panel), Ne-Ne (middle panel) and Ar-Ar (lower 
panel). 
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